* Reset settings and initialize log file
launch, path("share/hhstatus")

*-------------------------------------------------------------------------------
* Price and Wasserman (2024), "The Summer Drop in Female Employment"
*
* Description: Show heterogeneity in the summer drop by household structure.
*-------------------------------------------------------------------------------


* Prepare data and run models
*-------------------------------------------------------------------------------

if "$estimate" != "0" {
	* Loop over household structure (parental status, age of youngest child)
	foreach g in "h" "k" {
		if "`g'" == "h" {
			local gvar "hhstatus"
			local grange "1/4"
		}
		else if "`g'" == "k" {
			local gvar "youngest"
			local grange "0/3"
		}

		* Compute aggregate stocks and weights
		gzuse "$basepath/data/derived/cps_bms_sample.dta.gz", clear
		gcollapse (mean) emp (rawsum) wtfinl (first) month tmspline* weeks [pw = wtfinl], by(female `gvar' tm)
		tempfile core
		save `core'

		* Estimate specifications
		foreach f of numlist 0 1 {
			foreach `g' of numlist `grange' {
				* Reload sample
				use if female == `f' & `gvar' == ``g'' using `core', clear
				quietly tsset tm

				* Run specifications
				quietly ivreg2 emp ib5.month tmspline* weeks [aw = wtfinl], bw($bandwidth) robust small
				process_estimates, path("hhstatus") model("f`f'_`g'``g''")
			}
		}
	}
}


* Prepare estimates
*-------------------------------------------------------------------------------

* Load estimates into memory
load_estimates, path("hhstatus")

* Prepare coefficient labels for each month
make_coeflabels


* Plot seasonality by sex x household structure
*-------------------------------------------------------------------------------

foreach f of numlist 0 1 {
	if `f' == 0 {
		local flbl "men"
		local ymin "-3.5"
		local ymax "1.5"
		local ylabel "-3.5(0.5)1.5"

		local col1
		local col2
		local col3
		local col4
	}
	if `f' == 1 {
		local flbl "women"
		local ymin "-2.75"
		local ymax "1.25"
		local ylabel "-2.5(0.5)1.0"

		local col1 "color(black) ciopts(color(black) lwidth(vthin))"
		local col2 "color(gs4) ciopts(color(gs3) lwidth(vthin))"
		local col3 "color(gs8) ciopts(color(gs6) lwidth(vthin))"
		local col4 "color(gs8) mfcolor(white) ciopts(color(gs9) lwidth(vthin))"
	}

	local shadecolor "$ltgs"
	local blackcolor "black"
	local recast "recast(connected)"

	* Prepare background shading
	clear
	set obs 14
	gen rlow = `ymin'
	gen rupp = `ymax'
	gen rval = _n - 0.5

	* Left panel: seasonality by marital status
	#delimit ;
	coefplot
		(f`f'_h1, offset(-.15) recast(connected) `col1')
		(f`f'_h2, offset(-.05) recast(connected) `col2')
		(f`f'_h3, offset(+.05) recast(connected) `col3')
		(f`f'_h4, offset(+.15) recast(connected) `col4'),
		coeflabels(`coeflabels')
		title("{bf:By marital and parental status}")
		xtitle("")
		xlabel(, alternate)
		ytitle("")
		yscale(range(`ymin' `ymax'))
		ylabel(`ylabel', format(%5.1f))
		addplot(
			rarea rupp rlow rval if inrange(rval, 0.5, 5.5), color($ltgs) plotregion(margin(t=0 b=0)) below ||
			scatteri 0 0.5 0 13.5, recast(line) color(black) below)
		vertical
		legend(rows(2) size(*1) bmargin(0 0 0 0) keygap(1.5) order(
			- "No spouse present:"
			4 "No child <18"
			6 "Parent"
			- "Spouse present:"
			8 "No child <18"
			10 "Parent"))
		rescale(100);
	#delimit cr

	tempfile byhhstatus
	graph save "`byhhstatus'"

	* Right panel: seasonality by age of youngest child
	#delimit ;
	coefplot
		(f`f'_k0, offset(-.15) recast(connected) `col1')
		(f`f'_k1, offset(-.05) recast(connected) `col2')
		(f`f'_k2, offset(+.05) recast(connected) `col3')
		(f`f'_k3, offset(+.15) recast(connected) `col4'),
		coeflabels(`coeflabels')
		title("{bf:By age of youngest child}")
		xtitle("")
		xlabel(, alternate)
		ytitle("")
		yscale(range(`yrange'))
		ylabel(`ylabel', format(%5.1f))
		yline(0)
		addplot(
			rarea rupp rlow rval if inrange(rval, 0.5, 5.5), color($ltgs) plotregion(margin(t=0 b=0)) below ||
			scatteri 0 0.5 0 13.5, recast(line) color(black) below)
		vertical
		legend(rows(2) size(*1) bmargin(0 0 0 2) keygap(1.5) order(
			4 "No child <18"
			6 "Age 0 to 5"
			8 "Age 6 to 12"
			10 "Age 13 to 17"))
		rescale(100);
	#delimit cr

	tempfile bykidage
	graph save "`bykidage'"

	* Create two-panel figure
	graph combine "`byhhstatus'" "`bykidage'", l1title("            Difference relative to May (p.p.)", size(*.8)) imargin(vsmall)
	nicepdf "$basepath/output/hhstatus_`flbl'.pdf", panels(2) indirect replace
}


* Report the magnitude of the summer drop for select groups
*-------------------------------------------------------------------------------

* Report the summer drop for married women with children
estimates replay f1_h4

* Report the summer drop for women with a youngest child aged 6-12
estimates replay f1_k2

* Report the summer drop for women with a youngest child aged 13-17
estimates replay f1_k3


* Close the log file
unlaunch
